%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% exercise of CS
% sensing matrix is random 1 and -1.
% Although the number of -1 is small, it also can complete the recovery.

clear all; clc; close all;

N  = 128;    %% # of points
M  = 64;

nz = 16;

% generate Fourier basis F
F  = zeros(N,N);
for i = 1 : N
    for j = 1 : N
        F(i,j) = sqrt(1/N)*exp(-1j*2*pi/N*((i-1)*(j-1)));       
    end
end
phi  = F';

% X is the Fourier coefficent of signal f, it is sparse in Fourier domain F
X = zeros(N,1);
X(10) = 10; X(11) =  12; X(100) = 2;
% X(50) = 3;
% X(55) = -4;
% X(60) = -17;
% X(70) = 9;
% X(30) = 1;
% X(23) = 11;

% f is orignal signal
f = phi * X;


% Generate random 1 -1 matrix
psi = zeros(M,N);
idx = randperm(M * N);
psi(idx(1:nz)) = 1;

y  = psi * f;
A = psi * phi;


%[xhat, ~] = omp(A, y, 3);
 torr = 1e-3;
 cvx_begin
     variable xhat(N);
     minimize( norm(xhat,1) );
     subject to 
         norm(A * xhat - y, 2)< torr;
 cvx_end

figure;
hold on;
plot(X,'gx');
plot(xhat,'r+');
plot(f,'b.');
plot(y, 'r:');



%观测矩阵是A。程序中的意思是，y是待恢复信号，它在FH上是稀疏的（只有3个非0值）
%FH也叫表示基，representation basis
